Driven lattice gas of dimers coupled to a bulk reservoir 
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We investigate the non-equilibrium steady state of a one-dimensional (ID) lattice gas of dimers. 
The dynamics is described by a totally asymmetric exclusion process (TASEP) supplemented by 
attachment and detachment processes, mimicking chemical equilibrium of the ID driven transport 
with the bulk reservoir. The steady-state phase diagram, current and density profiles are calculated 
using both a refined mean-field theory and extensive stochastic simulations. As a consequence of 
the on-off kinetics, a new phase coexistence region arises intervening between low and high density 
phases such that the discontinuous transition line of the TASEP splits into two continuous ones. The 
results of the mean-field theory and simulations are found to coincide. We show that the physical 
picture obtained in the corresponding lattice gas model with monomers is robust, in the sense that 
the phase diagram changes quantitatively, but the topology remains unaltered. The mechanism for 
phase separation is identified as generic for a wide class of driven ID lattice gases. 



I. INTRODUCTION 

Driven lattice gas models have recently received much 
attention due to their possible application to intracellular 
transport lj. The connection of non-equilibrium trans- 
port to biology has been recognized almost 40 years ago 
by modelling ribosome motion on mRNA. In two pio- 
neering papers MacDonald et al. 0, Q have introduced 
a one-dimensional (ID) driven lattice gas model which 
later became known as the Totally Asymmetric Simple 
Exclusion Process (TASEP). 

In this model each particle occupies a site on a 
one-dimensional lattice and advances stochastically uni- 
directionally (hence the term "totally asymmetric"). The 
non-trivial properties arise by imposing the constraint 
that each site may at most be occupied by one particle, 
and that moves are forbidden if the target site is already 
occupied (hence the term "exclusion" ) . Of particular in- 
terest are systems with open boundary conditions where 
particles may enter and leave the ID lattice at the left or 
right boundary with rates a and (3, respectively. 

Remarkably, this ID driven system with open bound- 
aries exhibits a non-trivial phase diagram with a and /3 
as control parameters |4|. 

For low extraction rate f3 a high density (HD) phase 
occurs, where the density profile is controlled by the right 
boundary. Symmetrically, for low injection rate a a low 
density (LD) phase appears, where the density profile is 
determined by the left boundary. For high values of both 
a and [3 the bulk density is independent of the boundary 
conditions and the system carries its maximal current 
(MC). 

Exact solutions were found both for periodic and open 
boundary conditions through exact methods, such as re- 
cursive relations, Bethe and matrix product ansatz 0,0. 



The stochastic fluctuations have been investigated in 
terms of static and dynamic correlation functions by ex- 
act methods 0, 0, E} and phenomenological approaches 
[ToL lr"r| . The rich phenomenology, as well as the avail- 
ability of exact results established TASEP as on of the 
paradigmatic models in non-equilibrium statistical me- 
chanics [TllrllTii. 

Over many years the TASEP has served as a testing 
ground for a mathematical analysis of non-equilibrium 
problems. More recently, the application to biology reen- 
tered the focus of attention in the context of systems 
of extended particles. For example, ribosomes bind to 
mRNA occupying twelve codons and progressively ad- 
vance by a single codon through two cycle of GTP hy- 
drolysis (elongation step). The implications of mutual 
steric hindrance of £-mers (monomers, dimers, 12-mers) 
have been studied by stochastic simulations and 'refined 
mean-field theory' for ^-TASEP jHE^. The picture al- 
ready found for monomers has been corroborated also by 
exact methods (Bethe ansatz and mapping to zero 
range process [l8|). 

Intracellular transport, where a molecular motor (ki- 
nesin) is moving along a molecular track (microtubule) is 
yet another example of ID non-equilibrium transport in a 
biological system (see e.g. [ljj,|2(|). In marked difference 
to directed motion of ribosomes on mRNA here there 
is exchange of molecular motors with the surrounding 
cytosol acting as a particle reservoir. This observation 
has motivated a recently proposed extension of TASEP, 
by supplementing the unidirectional hopping with Lang- 
muir (on-off) kinetics |2lL l22j|. Single motor experiments 
suggest that both processes compete on the scale of the 
microtubule, i.e. a motor explores a significant fraction 
of the track before detaching j22| ■ To capture such an in- 
terplay mathematically, a mesoscopic limit has been sug- 
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gested, where local adsorption-desorption rates have been 
rescaled in the limit of large but finite systems, such that 
the gross rates are comparable to the injection-extraction 
rates at the boundaries The study of this system 

has shown an unexpectedly rich and new phase diagram 
with phase coexistence regions. Surprisingly in the meso- 
scopic limit the density profiles show a sharp discontinu- 
ity (shock separating low from high density phase) or 
cusp in a large portion of the phase diagram. 

Here we propose a model for intracellular transport 
incorporating both coupling to a reservoir as well as the 
finite extension of particles. The model is motivated by 
the fact that many processive molecular motors (e.g. ki- 
nesins, dyneins and myosin V) are composed of two heads 
(trail and lead) that bind specifically each to a subunit 
of the molecular track (a- and /3-tubulin in case of mi- 
crotubules). Hence, we concentrate on particles which 
occupy two sites, i.e. dimers, and advance by a single 
site per hopping event. 

The theoretical challenge arises from finding an appro- 
priate and quantitative description that captures both 
the on-off kinetics as well as the £-TASEP features on an 
equal footing. Whereas for monomers a straightforward 
decoupling of correlations yields a consistent mean-field 
theory, such an approach fails even to describe the chem- 
ical equilibrium of a pure on-off kinetics |l5j . In the fol- 
lowing paragraphs we review some aspects of our model 
where the geometric constraints usually play an impor- 
tant role introducing non-trivial static and dynamic cor- 
relations. 

At the static level the geometrical frustration emerges 
already in the equilibrium thermodynamics of the one- 
dimensional Tonks gas. Identical extended particles are 
distributed according to the grand canonical ensemble on 
a one dimensional lattice. A hard core repulsion reduces 
the task to a combinatorial problem |24l l25j . It turns 
out that the single site occupation densities do not de- 
termine the configurational probabilities, contrary to the 
monomer case. There the Langmuir kinetics specified 
in terms of attachment and detachment rates loa and 
ljd, or equivalently the binding constant K = ua/ojd- 
The coverage density is simply determined by a single 
site consideration, resulting in the Langmuir isotherm: 
pi = K/(K + 1) = lua/{ua + ^d)- This is to be con- 
trasted with the dimer case where a full evaluation of the 
grand canonical partition function is required leading to 
an equilibrium coverage density |2q: 



At a dynamical level a non-trivial approach to the 
steady state is found in the on-off kinetics in the fast 
attachment limit: at short time scales only deposition 
processes are frequent, while detachment processes are 
still unlikely. The kinetics belong to the class of prob- 
lems referred to as random sequential adsorption (RSA) 
and the so-called "Flory plateau" [2]} is reached expo- 
nentially fast. After this transient the particles can de- 



tach (on longer time scales) freeing new gaps; this pro- 
cess results in a long tail relaxation where the isotherm 
is approached with a power law i -1 / 2 . Such behavior has 
been explained using the mapping of the gap dynamics 
to a ID reaction diffusion system A + A — ► (for a re- 
view see [Isl H^, |23| and references therein) for which this 
anomalous power law dynamics is known. This two step 
relaxation has been found also in other models like in ir- 
reversible deposition of dimers (i.e. RSA) with diffusion 
[8ll | and is intrinsically related to the extended nature of 
the particles. 

Driving a system of extended particles far from equilib- 
rium, e.g. introducing the ^-TASEP dynamics, suggests 
a variety of competing time scales and ordering phenom- 
ena, i.e. non-equilibrium phase transitions. The purpose 
of this paper is to shed light on density and current pro- 
files in the stationary state for a model incorporating 
such aspects, as well as to provide a complete analytical 
theory supported by simulated data. 

This paper is organized as follows: in Sec. [H] we spec- 
ify the model under consideration and introduce the used 
formalism. Section ITTT1 exemplifies results of the stochas- 
tic simulations and puts extracted features into proper 
perspective. In Sec. IIVI we construct a suitable set of 
mean-field equations relying on results for £-TASEP and 
on-off kinetics. In Sec. with the help of these equa- 
tions, we rationalize the MC data and study the complete 
phase diagram. In the Conclusion, Sec. lVIl we discuss the 
robustness of the model, give a possible generalization to 
£-mers and summarize our main results. 



II. MODEL AND NOTATION 

We consider a finite one dimensional lattice with sites 
labeled i — 1, . . . , N, of unit length L — 1; consequently 
the lattice spacing reads a = L/N. The first two sites 
(i = 1,2) and the last two sites, (i = N — 1,N) rep- 
resent, respectively, the left and right boundary. The 
lattice is partially covered by dimers, i.e. composite ob- 
jects consisting of two monomers rigidly tied together. A 
dimer occupies two lattice sites, see Fig.^ we refer to the 
monomer to the right of the dimer as the lead head and 
correspondingly to the left monomer as the trail head. 

A microscopic state C of the system consists of a con- 
figuration of non-overlapping indistinguishable dimers on 
the lattice. Non-overlapping requires that no two heads 
(trail or lead) occupy the same lattice site. We specify 
the dynamical evolution by the following set of updating 
rules: 

• If the lead head of a dimer occupies a site i = 
2, . . . ,N — 1 and the following site is empty, the 
dimer advances one step with unit rate. 

• At site i — 1 a dimer enters the lattice with its trail 
head with rate a provided that the first two sites 
are empty. 
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• A dimer with its lead head on site i = N leaves the 
lattice with rate (3 emptying the last two sites. 

• Everywhere in the bulk (i.e. trail head on a site 
i = 2, . . . , N—2) a dimer leaves the lattice emptying 
two sites with a site independent detachment rate 

LU D . 

• Everywhere in the bulk (i = 2, . . . , N — 2) a dimer 
enters the lattice with its trail head, provided that 
the considered site as well as its right neighbor 
is empty, with a site-independent attachment rate 

LU A . 

The first three rules encode the usual Totally Asym- 
metric Simple Exclusion Process of extended objects (£- 
TASEP), while the last two implement the coupling of 
the lattice to a reservoir of dimers with fixed chemical 
potential |3^ |. 

The rules entail that a complete description of the un- 
derlying stochastic Markov process is given in terms of 
the time-dependent probability P(C,t). The time evolu- 
tion is governed by the associated master equation 

d t P{C) = [W c ^cP(C, t) - W c ^c'P(C, *)] , (2) 

where the transition rates Wc—>c can be inferred form 
the dynamical rules. 

It is expected that the system evolves towards a sta- 
tionary ensemble P st (C) for long times. Due to the lack 
of detailed balance this macroscopically stationary state 
does not correspond to an equilibrium ensemble, i.e. there 
is no Gibbs-Boltzmann measure to be inferred from ther- 
modynamics arguments. The macroscopic quantities we 
aim to compute are understood as averages over this dis- 
tribution. 

To exploit probabilistic methods we label sites accord- 
ing to: empty (state s = 0), occupied with the trail head 
(state s = 1) or with the lead head (state s = 2). The 
configuration C can be represented as a string of occu- 
pation number C = {rii, . . . , tin}, with m — {0,1,2}. 
The quantity of fundamental interest in the following is 
the average site-dependent dimer density in the station- 
ary state, In particular we shall derive equations for the 
density of the lead head (m — 2) at site i, pi = {8 ni ^)- 
Since the dimers are rigidly connected, the corresponding 
density for the trail head at site i simply reads Pi+x- 

The probability of having site i in state s will be de- 
noted by p(i, s) = P(rii = s). Since the states mutually 
exclude, the probability p(i, 2) immediately yields the av- 
erage lead head density: pi = p(i, 2). The coverage den- 
sity pi is sum of the lead head density and the trail head 
density, i.e. pi — p(i, 2)+p(i+l, 2). Consistently with the 
notation the joint probability will read p(i, s; j, s'). We 
also shall need the conditional probability p(i + 1, s\i, s'), 
i.e. the probability for site i + 1 to be in state s provided 
site % is in state s' . 

In the simulations we will concentrate on the average 
coverage density p°, and the spatially resolved current ji, 



defined as the flux of particle through site i per unit time. 
Only averages on the stationary state will be considered: 
a comparison between ensemble and time-moving aver- 
ages corroborates the hypothesis that the system is er- 
godic. 
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FIG. 1: Schematic representation of the model and allowed 
moves: forward jump (with rate r = 1), entrance at the left 
boundary (with rate a), exit at the right boundary (with rate 
13), attachment (with rate to a), and detachment (with rate 
lod) in the bulk. 



As indicated in the Introduction, a non-trivial compe- 
tition between boundary induced phenomena and bulk 
dynamics arises in the mesoscopic limit, i.e. the gross 
on/off rates, SIa and fin, should be of the same order 
of magnitude as the boundary entry/exit rates a and (3. 
This requires that the single site on/off rates, uja, <^d, 
scale with the system size N as: 



UJA 



N ' 



N 



(3) 



with Qa an d &d fixed. 



III. MC RESULTS AND DESCRIPTION OF 
THE PHASE DIAGRAM IN GENERAL 

In this section we present the method used to per- 
form numerical simulations and the results obtained from 
them. Density profiles much different from the TASEP 
motivate our interest in the model. 



A. Method 

We have simulated the stochastic equation through a 
kinetic Monte Carlo algorithm, to determine the average 
density profile in the stationary state with high accuracy. 
The results have been obtained using a random sequen- 
tial updating algorithm by Bortz, Kalos and Lebowitz 
(BKL-method or n-fold way) |^ : a list of all sites 
which are possible candidates for a successful move is 
stored and dynamically updated. The method is (for the 
present case) faster than the conventional algorithms and 
constitutes a reliable way to simulate real time dynam- 
ics, although here we focus on equal-time averages in the 
stationary state (note that the algorithm is equivalent to 
Gillespie's, commonly used in chemical reactions 

Explicitly, we keep three lists of sites for the respec- 
tive allowed moves: there are Nj sites where particles 
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can perform a jump forward, Na sites that may accept 
a particle from the reservoir, and Nn sites from which 
a particle can detach. Furthermore we keep track of the 
occupation of the boundaries (no and njv). One of these 
moves is selected randomly with the appropriate weights: 
lxJVj, ujaxNa, LUDxNo,ax (1 — ni), /3xnjv-i- Then we 
increment time by an interval drawn from an exponential 
distribution with timescale T = (a(l — no) + /3njv + Nj + 
u>aNa + ujdNd)^ 1 - The lists are updated correspond- 
ingly and the procedure is iterated for many events. We 
have started the algorithm from an empty configuration 
and after the stationary state has been reached, the cov- 
erage density and the current have been measured. Mov- 
ing time averages typically cover a window of O(10 6 ) time 
units. We have investigated finite-size effects by consid- 
ering lattices varying from 128 up to 4096 sites. Fur- 
thermore, we have checked the ergodicity by comparing 
time-moving and ensemble averages. The longest simu- 
lations have been performed approximately a day on a 
1GHz processor. 



B. Results 

We exemplify our results by fixing the binding con- 
stant to K = 3 and the entrance/exit rates to a = 0.1 
and f3 = 0.6; see Fig- EI We have measured the stationary 
density profiles as well as the current for different values 
of the detachment rate f2 d ■ For very small 0^ we find al- 
most flat profiles for both p and J, in accordance with the 
picture of simple dimer TASEP [lH[l]|. On the other- 
hand, for very large detachment rates the on/off kinetics 
dominates the bulk of the density profile. However, con- 
trary to equilibrium thermodynamics of on/off kinetics, 
a large directed current is simultaneously present. 

For intermediate values of the detachment rate, neither 
the density profile nor the current is constant anymore. 
There appear rather narrow regions where the density 
steeply rises from an almost linear profile to one that 
approaches the isotherm value. Such findings are similar 
to the ones of monomeric TASEP coupled to Langmuir 
kinetics |22l l3^ | , where it has been shown that the steep 
increase corresponds to a domain wall of a coexistence 
phase in the mesoscopic limit. 



IV. CONSTRUCTION OF THE MF EQUATION 

Since the stationary density profiles for TASEP of 
monomers coupled to Langmuir kinetics have been suc- 
cessfully described in terms of mean-field theory, it ap- 
pears promising to look for an analytical description for 
dimers also. However, it turns out that in the case of 
dimers the conventional mean-field approach, i.e. a de- 
coupling of correlations at lowest order starting from the 
quantum Hamiltonian formulation, see e.g. |39j . is inap- 
propriate: the non-trivial correlations, arising from the 



p c 0.4 







I n D = i 

. ! 


/ / Or, =0.1 1 
y Q.D = 0J)5^~ 


Q.D = 0.001 



0.2 0.4 rj. 0.6 0.8 1 



0.05 



n D = 10 

/ Q D = 1 


= 0.2 






f/y 


/VLd = 0.1 




= 0.05 








= 0.01 








= 0.001 



0.2 0.4 X 0.6 0.8 1 



FIG. 2: Stochastic simulations results for stationary density 
profiles and corresponding currents for a system of N — 4096 
lattice sites. Parameters have been fixed to a — 0.1, (3 = 0.6, 
K — 3, and Q.d are varied as indicated in the figures. 



fact that the particles are extended, manifest themselves 
already at the level of pure TASEP of dimers . 

The failure of conventional mean-field approach re- 
quires the introduction of phenomenological methods. 
The aim of this section is to identify useful strategies that 
are capable to deal first with the simple on/off kinetics 
of extended objects without directed motion and second 
with pure TASEP for dimers. We will then combine the 
two cases to obtain a formulation for the complete prob- 
lem, and we will provide arguments that the aspects of 
correlations are captured on a reliable level for both pro- 
cesses. 

The density profile in the mesoscopic limit will be de- 
scribed in terms of a balance equation. The evolution of 
the course-grained dimer density is determined by: (a) a 
flux j of dimers within the lattice that encodes the asym- 
metric diffusion, (b) two source terms that represent the 
on-off kinetics of dimers: 

dtp = -d x j + F A -F D . (4) 

The non-trivial task is to provide reasonable expressions 
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for the current and the sources in terms of the local den- 
sity. In Sec. |V] we shall compare the analytical results 
to stochastic simulations and demonstrate the validity of 
this approach. 



kinetic rates Qa and ftp we obtain the rate equation in 
the mesoscopic limit: 

d T p = fl A (1 ~ 2p) - tt D p + 0{a 2 ) . (9) 

1-/3 



ON-OFF kinetics of dimers 



1. Probabilistic approach 



The dynamics of the on-off kinetics in a rate equation 
approach is described by the difference between a gain 
term F, A and a loss term FP: 



It is natural to expect a density limited to the range 
[0, 1/2] , because there cannot be more than one lead head 
per site. In the stationary regime the physical solution 
corresponds to the constant density pi, Eq. QJ, deter- 
mined solely by the binding constant K = CIa/CId- This 
value of the isotherm is consistent with the one obtained 
by McGhee and Von Hippel f26| for the general case of 
£-mers. 



d tPl =F t A -Fl 



(5) 



These two terms are proportional to the attachment and 
detachment rates, loa and wp. Clearly the detachment 
flux is exactly given by Ff = a; n p.- . |36| 

An attachment event at site i (lead head) requires 
both sites i and i — 1 to be empty (i.e. in state (0)). 
The corresponding attachment flux then reads F$ = 
uJAp{i — 1,0; i,0), and knowledge of the joint probabil- 
ity is needed. In order to obtain closed equations for the 
lead head density we break the hierarchy of joint proba- 
bilities resulting from the master equation. Here, we rely 
on an ad hoc approximation for the conditional probabil- 
ities: p(i,0\i-l,0) ~p(i,0|(i-l,0)v(i-l,2)).|33 The 
idea is that the probability of having site i empty is inde- 
pendent of having a lead head or a hole in the site i — 1: 
note that this is precisely what usually the mean-field ap- 
proximation does in monomeric lattice gases, where cor- 
relations are broken factorizing the density: (71^+1) = 
(rii) (rii + i) . In Appendix lAl we show that this assumption 
yields a closed expression for the joint probability: 



p(*,0;i-l,0) 



p(i-l,0)p(i,0) 
l-p{i,2) 



(6) 



Eliminating the probabilities for single empty sites and 
collecting results, the rate equation for the lead head den- 
sity reads: 



dt P i = LOA 



(1 - pi-i - pi)(l 



Pi+i) 



1 



UJDPi 



(7) 



As a next step we perform a continuum limit. Rela- 
belling position as fraction of the system size, x = i/N, 
the density pi becomes a function p(x) in the mesoscopic 
limit N — ► 00. The average density at the neighboring 
sites are obtained by expansion in powers of the lattice 
constant a — 1/N: 

p(x ± a) = p(x) ± ad xP {x) + -a 2 dlp[x) + 0{a 3 ) . (8) 

Introducing a new rescaled time t — at and eliminating 
the single site on/off rates loa and lod in favor of the gross 



2. Tonks gas approach 

To corroborate the results presented above, we com- 
pare the isotherm of the rate equation with the one com- 
puted from the grand canonical partition function of the 
Tonks gas. In particular, we shall identify the fugacity z 
with the binding constant K. 

The statistical mechanics of the Tonks gas on an open 
lattice consists of distributing extended objects on the 
N sites, respecting the excluded volume constraint. In 
the case of the canonical ensemble this reduces to the 
combinatorial problem of counting the number of ways to 
distribute n dimers on N lattice sites. Using the standard 
trick to represent the two occupied sites of each dimer by 
a "stick" and the N — 2n empty sites by "balls" [2{| , one 
immediately concludes 



Z(n, N) 



N - 



(10) 



The corresponding grand canonical partition function is 
readily evaluated using the previous result 



N 



Z N = ^2z n Z(n,N) 



n=0 



( v T+4i+i) Ar+1 + (VT+4i-i) A,+1 

2 Ar + 1 v / TT4i 



(11) 



where z denotes the fugacity 38] . The average density is 
then obtained in the thermodynamic limit by 



1 ( 1 

p = hm — z—— In Zm = — 1 ■ 

P n^oo N dz 2 V V4ITT 



(12) 



Since the master equation of pure on/off kinetics ful- 
fills detailed balance, the stationary state is given in 
terms of the grand canonical ensemble. By a single de- 
tachment event a configuration C of n dimers connects 
to some new configurations C' of n — 1 dimers. For 
such configurations the detailed balance condition im- 
plies P(C) / P(C) — loa/lod = K. Conversely, the ratio 
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can be determined also from the grand canonical Boltz- 
mann factors. Since, apart from the total exclusion, en- 
ergy does not enter the problem, the probabilities are de- 
termined by the number of dimers only: P(C) = z n / Zjq 
and P(C) = z n ~ 1 /Zw- Combining both expressions, we 
conclude z — K as in the case of monomers. 

The last result shows that the "isotherm" obtained us- 
ing ad hoc approximations fsub-Section llV A lfl is at least 
consistent with thermodynamics. Let us mention that in 
equilibrium the situation is probably better than in the 
general case, since correlations are quickly washed out 
by the coupling to the reservoir. The relaxation towards 
equilibrium can be rather different from the naive picture 
of rate equations, see e.g. |4Cj . 



B. TASEP of dimers 



p(l) = ^. (17b) 

The details are presented in Appendix [5] The prop- 
erties of these equations have been studied in details 
in Ref. ^||- The resulting phase diagram is topolog- 
ically equivalent to the one of TASEP of monomers. 
The continuum analogue for the current-density relation, 
Eq. lfT5|). reads 

1 - P 

and exhibits a maximum at a distinguished density p* = 
1/(2 + \/2) ■ This value will play an important role in the 
more general case of TASEP coupled to on/off kinetics. 



In pure TASEP particles cannot leave or enter the 
track except at the boundaries. Correspondingly in the 
bulk a conservation law holds 



dtP i = ji-i ~ ji 



(13) 



The currents ji can be determined by probabilistic argu- 
ments. A hopping event of a dimer at site i requires first 
that site i is occupied by a lead head and second that site 
i + 1 is empty. Since we fixed the hopping rate to unity 
this yields 



ji=p(i,2;i +1,0). 



(14) 



Hence, also in the case of pure TASEP the master equa- 
tion induces a whole hierarchy of joint probabilities. To 
close the equation we rely on the same ad hoc approxi- 
mation for the conditional probabilities as for the on/off 
kinetics. In the appendix we show that this implies for 
the current a closed expression 



pi(l - p l+ i - pi +2 ) 



(15) 



1 - Pi+l 

Performing the continuum limit as in Eq. @ one obtains 



d T p = -d x 



p(\ - 2p) 



1 



1 - 2p 2 a 

— w^ xP 



(1 



0(a 2 ) 



.(16) 



Since, at the very end, we consider large but finite sys- 
tems we have kept the leading correction in a, which will 
turn out to be relevant for the formation of shocks and 
boundary layers. 

The equation of continuity in the bulk has to be sup- 
plemented with appropriate boundary conditions. Rely- 
ing again on the same ad hoc assumption for the condi- 
tional probability, and performing the continuum limit, 
one finds 



P{0) 



a 



1 + a 



(17a) 



C. Combining TASEP and on/off kinetics 

As anticipated at the beginning of this section, we iden- 
tify the current and the source terms of Eq. Q and we 
re-expressed everything in term of the density. The cur- 
rent is described by Eq. (|18fl and the on/off kinetics by 
Eq. @. By taking the mesoscopic limit we have intro- 
duced the gross attachment and detachment rates and 
rescaled the time consistently in both processes. The 
balance equation to order a can be easily rewritten as 



d T p = -d x 



a 

2(1 



d x p p{l - 2p) 



P? 



(l-2pf 

(1 - P) 



(19) 



with the boundary conditions, Eq. (|17fl . These equations 
constitute the mean-field description to be discussed in 
the following section. 

Equation (|19fl may also be recast in a form reminiscent 
of a non-linear Smoluchowsky equation |4l| : 



d T p = -d x [A(p)p - B(p)pd xP ] + C(p) 



(20) 



where the coefficients A, B and C depend on p. The 
second coefficient will prove not to be relevant since it 
scales to zero with increasing the system size to infinity. 

We may now interpret j(p) = pA(p) as an effective 
particle current and C(p) as a source/drain term; see 
Fig. 13 for an illustration. As will become clear in the 
following section the maximum of j(p) and he zeros of 
C(p) play an important role for the analytical from of 
the stationary solution. 



V. SOLUTION OF THE MEAN-FIELD 
EQUATION AND PHASE DIAGRAM 

This section is devoted to a discussion of the phase 
diagram emerging from combining TASEP of dimers with 
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FIG. 3: (a) Graph of the current-density relation pA(p) in 
Eq. 1201 . (b) Graph of the source-drain C(p) term in the 
Smoluchowsky equation (1201 for different values of K and 

n D = o.i. 



on/off kinetics. We shall derive density profiles within a 
mean-field approach, identify different phases and phase 
boundaries as the kinetic rates are varied. In particular, 
we demonstrate the role of the binding constant in the 
topological changes of the phase diagram as a certain 
critical value is reached. 



A. Emergence of shocks and boundary layers 



shock one obtains 
j(x w + Sx) - j(x w - Sx) 



x w -\-5x 



x w —5x 



(F A -F D )dx. (22) 



In the limit a \ the l.h.s. reduces to j a (x w + Sx) — 
j/3(x w — Sx) (where we have defined j a — jO Q (l — 2p a )/(l — 
p a ) the current set by the left boundary and similarly for 
the right one jp). Performing the limit Sx \ the r.h.s. 
of Eq. I|22|l vanishes and yields the matching rule in terms 
of the currents 



(23) 



Equivalently the rule implies for the densities at the 
matching point 



Pa i.^w) 



1 1 - 2p p (x w ) 

2 l-pp(x w ) 



(24) 



Let us make a comment. The fact that the in-track cur- 
rent is continuous at the shock is consistent with the idea 
of the mesoscopic limit. The fluxes due to attachment 
and detachment are important only on the length scale 
of the system size. Locally, i.e. on the scale of the lattice 
constant, the balance equation is dominated entirely by 
the unidirectional hopping process (TASEP). 



The mean-field equation for the stationary density pro- 
file, Eq. H19|) supplemented with d T p = 0, is a second 
order differential equation as long as the lattice constant 
a is small but finite. In the limit a \ the differential 
equation degenerates into a first order one. However, the 
density profile has to match the two boundaries, and it 
appears that the problem is overdetermined. Indeed one 
can construct separately a density profile fulfilling either 
the left or right boundary condition. The two branches 
will be referred to as p a and pp, respectively. In general, 
the two branches do not join smoothly together. One ex- 
pects the solution of the second order differential equa- 
tion to be well approximated by the two branches p a and 
pp except for a region where the solution rapidly crosses 
over from one branch to the other. Upon decreasing the 
lattice constant a the crossover region shrinks leading 
eventually to a discontinuity of the density profile for 
a = 0. Depending on the position of such a discontinuity 
x w we refer to it as a shock (0 < x w < 1) or a boundary 
layer (x w — or x w = 1). 

To locate the position of the shock or boundary layer 
we rely on a very general argument already used in 
[22I l39j . The second order differential equation for the 
stationary density profile can be written in the form of 
a balance equation d x j — Fa — Fd, where the current 
along the track reads 



a d x p p(l - 2p) 



(21) 



2(l-p)a 1-p ■ 
Integrating over a small region of width 25x around the 



B. Analytical solution 

The left and right branches of the stationary density 
profile are determined by solving Eq. (|19fl . once we set 
the l.h.s. to zero and discard the second order derivative, 



-d x 



p(l - 2p) 



1-P 



(1_2^ 



(25) 



(where the binding constant K = f^/f^ has been intro- 
duced) which has to be supplemented by the appropriate 
boundary conditions Eq. 117fl . By separation of variables 
the general solution G{p) = x + const is obtained after a 
straightforward integration with 



G{ P ) 



1 



p)+A + \n(p+ -p) 



+A_ln|p_-p|], (26) 
where we introduced the zeros of the source term 



1 ( 1 
P± = - I 1 ± 



and the amplitudes 



-4, 



K 



VTTIk 



1 - 4p ± + 2p\ 



l + 4K(l-p ± )(p ± -p^) 



(27a) 



(27b) 



The function G(p) exhibits a singularity in the physical 
regime 0</9<l/2at the isotherm of the on/off kinet- 
ics pi = p-. One easily checks that G(p) exhibits an 
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extremum when the corresponding current, Eq. (|18fl . is 
maximal, i.e. at p = p* . For pj > p* (pj < p*) the ex- 
tremum of G(p) corresponds to a maximum (minimum). 
A change of topology of the graph of the solution oc- 
curs once both distinguished densities coincide pj = p*, 
which happens at a critical value of the binding constant 
K* = (1 + \/2)/2. At this special value the amplitude 
A_ vanishes and the function G(p) becomes smooth and 
monotonic in the physical regime < p < 1/2. 

After matching the boundary conditions, the left and 
right solutions are obtained up to inversion of a function 



L a (p a ) = G(p)-G(a/(l + a))=x, (28a) 



R^pp) ee G(p) ~ G((l - P)/2) + l=x. (28b) 

Upon inverting, the singularity of G(p) transforms into a 
horizontal asymptote, whereas the extremum manifests 
itself as a branch point, see Fig. 0J The functions are 
multivalued with three branches, to be referred as Wq, 
W ~andW_i. 



C. High binding affinity: K > K* 

For large binding constants, K > K* , the equilibrium 
density of the on-off isotherm is larger than the critical 
density at which the current of pure TASEP is maximal, 
pi > p*, see Fig.|U 

For densities smaller than pi, the track accumulates 
particles from the reservoir, since 



P 0.5 



Fa Fj 



D 



l + AK 



D~ 



1 



[{P ~ P-){P ~ P+)\ 



(29) 



Hence the isotherm acts as an attractor. This property 
has to be reflected in the left and right solution solving 
the balance equation for the current, d x j = (d p j)(d x p) = 
F A - F D . For p(0) = a/(l + a) > p* , the branches 
are unstable. For densities p > pi the coupling to the 
reservoir induces a decrease in density, whereas Wq is 
increasing. Similarly for p* < p < pi the net gain due 
to the on/off kinetics is not reflected by the decreasing 
function Wq . One concludes that for the left solution the 
only physical branch is W-%. For a > a* = p*/(l — p*) = 
1/(1 + ^2), i.e. the critical entrance rate of dimer TASEP, 
a left boundary layer necessarily appears. 

For the right solution the argument is just reversed, i.e. 
the physical branches are Wj^ . The critical exit rate f3* 
of dimer TASEP is determined by p(l) = (l-/3*)/2 = p* 
and equals (3* = 1/(1 + V2). For exit rates larger than 
the critical one a right boundary layer emerges. The 
complete phase diagram is obtained by identifying dif- 
ferent scenarios to match left and right solutions. The 
boundary conditions, Eq. i|17|l can be satisfied by the 
stable branches of the left and right solution only for 




0.4 x 0.6 



P 0.5 




0.4 x 0.6 



FIG. 4: The functions W-x (dashed lines), W " (full lines), 
Wo + (full lines) for (a) K = 4, Q D = 0.05 (high binding 
affinity) and (b) K — 0.25, Qd = 0.05 (low binding affinity). 
The branches Wq and Wq are separated by the isotherm 
density pi . Varying a or /3 results in a simple shift parallel to 
the a;— axis. In Fig.(a/b) the thick lines represent the solutions 
where the branch point hits the right/left end of the system. 
The dotted lines indicate the critical density value p* (where 
branch point occurs) and the isotherm pi. 



a < a*, (3 < f3* respectively. The merging of the two 
solutions by Eq. (|24|l is achieved at some position x w . 
Depending on whether the (i) < x w < 1, (ii) x w < 
and (iii) x w > 1, a domain wall in the bulk emerges (i), a 
boundary layer at the left (ii) or right (iii) appears. For 
(3 > (3* an additional boundary enters the system. The 
result of such an analysis is presented in Fig. and the 
chain of arguments follows the one of Ref . [39j . 

For example, for large entrance rate a and low exit 
rates (3 the density profile is in a high density phase 
(HD), i.e. the solution p(x) exceeds p* , and is given by 
the stable right solution pp{x). At the left end x = 
a boundary layer appears. Similarly a region of a low 
density phase (LD) and low density-high density phase 
coexistence (LD/HD) with appropriate boundary layer 
and domain wall is inferred from the figure. 
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FIG. 5: High binding affinity; phase diagram for K = 2, 
Q.D = 0.1. In the low-density (LD) phase the profile is given 
by the left solution p a and a boundary layer appears at the 
right (r). The position of boundary layers (r, I), as well as 
the relevant bulk solution (p a , pp, and pp*) are indicated for 
all phases (HD: high density phase, LD/HD: low density-high 
density coexistence phase, M: Meissner phase, LD/M: low 
density-Meissner coexistence phase, see text). 



For (3 > (3* the right solution is replaced by the critical 
right solution pp=p* , which constitutes the analogue to 
the maximal current phase in TASEP : the system indeed 
cannot transport more than the current jp* imposed by 
the critical exit rate (3* . 

Changing the exit rate beyond (3* does not affect the 
density profile except in a small boundary layer. This 
phenomenon is similar to a superconductor where an ex- 
ternal magnetic field does not enter the sample except 
for a short penetration depth. This analogy suggests to 
denote this region in the phase diagram by M (Meiss- 
ner) phase. The coexistence of low-density and Meissner 
phase (LD/M) extends the coexistence phase (LD/HD) 
for /?>/?*. 

The phase boundaries can be computed analytically up 
to inversion of a function. For example, the HD - LD/HD 
phase boundary is obtained by requiring the domain wall 
to fall on the left end of the system x w = 0. Thus the 
continuity condition, Eq. (|24|l translates the boundary 
condition p(0) = a/(l + a) to p(x + 0) = (1 - a)/2. In 
terms of the implicit solution, Eq. J28bj, this implies 



by placing the domain wall at x v 
1 + /3 



= 1, which yields 

1. (31) 



For (3 > (3* the preceding equations still hold provided (3 
is replaced by the critical exit rate /3*.|43|] Correspond- 
ingly, phase boundaries degenerate to straight vertical 
lines. 

Density profiles exemplifying different regions of the 
phase diagram are showed in Fig. and compare nicely 
to the corresponding simulated data. 
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FIG. 6: (a) Density profiles for K = 2 and Q. D = 0.1 
in systems of N = 4096 lattice sites; parameters (a, (3) = 
(0.1, 0.6), (0.2, 0.6), (0.6, 0.6), (0.6, 0.2). (b) Approach of the 
mesoscopic limit of the density profile for a — f3 — 0.1, 
f2.D = 0.1, K — 2 and different system sizes (TV = 
128,256,512,1024,4096); the domain wall becomes steeper 
with the increasing of N. Wiggly lines represent simulated 
data, while dotted lines indicate analytic solutions 



Rp 0-^=0. (30) 
Similarly the LD- LD /HD phase boundary is determined 



The low density phase disappears from the phase dia- 
gram upon increasing flu at fixed K. Such a topological 
change occurs, once the phase boundary of the LD/M 
coexistence phase hits the /3-axis of the phase diagram 
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(i.e. a 
reads: 



0). The condition is inferred from Eq. (|31|l and 



L a =o 



1 + 13* 



1. 



(32) 



A numerical solution of the preceding equation for Q* D as 
a function of K, is presented in Fig. [7| 

In the limit FId — > + the LD/HD coexistence phase 
shrinks continuously to the line a = (3 and one recovers 
the dimer TASEP phase diagram. For Cl D > fl* D , the 
coexistence phases constitute only a marginal region in 
the a — (3 plane, located close to the /3-axis. The phase 
diagram is dominated by the HD and M phases, however 
the density profile approaches the constant value pi in 
the bulk, as expected from pure on-off kinetics. |44| 



matching the appropriate left and right solutions. A cut 
of the parameter space at fixed £Id and K is presented 
in Fig. [7| The high density phase disappears from the 
phase diagram for detachment rates exceeding a critical 
value Q* D (K), see Fig-EllH 

Let us stress that there is an important difference with 
respect to the monomeric analogue, namely the particle- 
hole symmetry is no longer present. Correspondingly, a 
symmetry transformation does not yield the phase di- 
agram of the "low affinity" from the one of the "high 
affinity" . Yet, the exchange of high densities and low 
densities, left and right, exit by entrance rates yields the 
correct topology of the phase diagram as well, as quali- 




FIG. 7: Critical detachment rate fl* D as a function of the 
binding constant K (solid line shows Q.* D (K) for dimers, dot- 
ted line for monomers). The dashed line K* separates the 
cases of "high" form the "low" binding affinity (for dimers). 
For large K the critical detachment rate approaches zero. In- 
sets represent the different topologies of the a — (3 section of 
the phase diagram. 
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FIG. 8: Low binding affinity; phase diagram for K = 0.8, 
Q.D = 0.1. Legend as in Fig. 



E. The critical case K = K* 



D. Low binding affinity: K < K* 

If K < K* the functional form of the solutions is still 
given by Eqs. (|26l28al28b|l . However, now the isotherm 
p = pi lies below the critical density, pi < p* . By the 
same reasoning as above, see Eq. (I29II . the left stable 
branches are the one which approach the isotherm, W^ 1 
in Fig.^J), while the only stable solution for the right side 
is W-i, i.e. reversed to K > K* . 

Using the same methods as for the high binding con- 
stant, phases and phase boundaries are identified by 



The limiting case when the binding constant equals the 
critical one, K — K* , requires a separate analysis, since 
the amplitude A_ vanishes, cancelling exactly the loga- 
rithmic divergence at p = p_ (nota bene p_ = p* = pi). 
This signals that p = p* = const is a spurious solution of 
the balance equation for the stationary density profile. It 
corresponds to the maximal current (MC) phase of pure 
dimer TASEP. Indeed for K = K* Eq. simplifies to 



(P-Pl 



2p_(2W2) / 1 



= 0. 



(33) 



11 



The non-trivial solution Eq. (126[l specializes to 



^^(inU-ri-^^-'))'' 34 ' 

which is a monotonic function in the physical regime. 
Correspondingly, the inverse function is single valued. 

The properties of the phase diagram for K = K* 
are obtained by assembling the resulting left and right 
branches, as well as the constant critical density p* , for 
different boundary conditions. We have determined nu- 
merically the line where the domain wall hits the left and 
right end, x w = and x w = 1, respectively by solving 
Eqs. i|30l and l^T|) . Moreover both left and right solutions 
now can approach the critical one (i.e. p = pi) in the 
bulk and therefore a triple phase coexistence LD-MC-HD 
is possible. This phase is bounded by the lines f3 = (3* 
above and a — a* on the right. Furthermore the curve 
g(a,P), where the left and right currents match addi- 
tionally the critical one, separates the triple phase coex- 
istence from the LD-HD phase. Implicitly g(a, j3) is given 
by the condition L a (/3*/(l + /?*)) = Rp ((1 - a*)/2). 

Depending on the value of the kinetic rate Qd, different 
topologies for the phase diagram arise. For example, for 
low rip seven phases are present in the a— (5 diagram (see 
Fig. Upon increasing JTq, a first topological change, 
analogous to the case of high binding affinity, occurs, 
i.e. the LD phase exits the phase diagram. By solving 
Eq. (|32[1 . one obtains the critical value fi* w 0.33945, see 
Fig. □ At another critical value ft£ w 0.52496 the HD 
phase disappears, and at still higher detachment rates 
flo > fi* = 0.86441 the LD-HD coexistence phase is no 
more present and the phase diagram (but not the density 
profiles) is independent of f2£>. Representative graphs of 
these topologies are shown in Figs. 171 and ITUl 



VI. CONCLUSIONS 

We have studied the non-equilibrium dynamics of a 
ID driven lattice gas with open boundaries, where indi- 
vidual particles are considered as dimers occupying two 
lattice sites simultaneously. Upon combining a general- 
ized mean-field theory with extended stochastic simula- 
tions, we have explored the phase behavior of the non- 
equilibrium steady state as a function of the kinetic rates. 
In particular, we have analyzed the interplay between the 
coupling of the system to a reservoir at the boundaries 
and at the bulk. 

In the following we compare our results with those 
obtained for the analogous problem of simple TASEP 
of monomers coupled to Langmuir kinetics. The bal- 
ance equation in the stationary state that constitutes the 
equivalent of Eq. (|T§)| has been derived in jSJ H^j : 



1 



d x [p(l-p)]=K(l-p)-p. 



(35) 
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Although the structure of the equation is the same, 
d x j = Fa — Fd, the current and the source terms obey 



FIG. 9: Critical case; phase diagram for the case K = K* , 
Q,d =0.1. Seven phases appear: LD, MC, HD and all the pos- 
sible coexistence: LD/MC, MC/HD, LD/HD, LD/MC/HD. 
Legend as in Fig. |5] 



a rather different functional dependence. In particular, 
the equation exhibits a particle-hole symmetry, i.e. den- 
sity profiles and the complete phase diagram is symmetric 
with respect to the following transformation rule p(x) ^ 
l-p(l-x), n D «-> A , a <h> [3. Also, a TASEP for dimers 
disregarding the on-off kinetics can be cast in a particle- 
hole symmetric form, introducing an effective hole den- 
sity p by the requirement j(p) = j(j>) and interchanging 
a with /3; this explain why the phase diagram is symmet- 
ric along the diagonal. Similarly, a pure on/off kinetics is 
invariant with respect to the transform Qa &d and an 
effective hole density p = (1 — 2p) 2 / (1 — p). However, one 
may easily show that the two symmetries of simple dimer 
TASEP and on/off kinetics never coincide for any K > 0, 
implying that the coupling of both processes breaks the 
particle-hole symmetry. 

It is remarkable that the topology of the resulting 
phase diagram is almost unaffected by this broken sym- 
metry. In fact, the single qualitative change of topology 
occurs at K = K* for intermediate tip, see Figs. [7| and 
110b . The origin of the robustness of the picture can be 
traced back to the following observations: (i) in a mean- 
field approach the balance equation for the stationary 
density profile is an autonomous first order ordinary dif- 
ferential equation; (ii) the current-density relation ex- 
hibits a single extremum for both cases [Ifj; (in) the so- 
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FIG. 10: Different topologies of the phase diagram for K = K* and (a) Sl^ < < fit, (b) < Qd < fic an d (c) Q,d > 
Phases are indicated in the graphs. 



lution p(x) is obtained by quadrature and inversion of a 
function resulting in a multivalued function with a simple 
branch point (due to the maximum in the current-density 
relation) and a horizontal asymptote (due to the on/off 
isotherm). Although the shape of the explicit solution 
depends on the details of the density current relation 
and the functional form of the on/off kinetics, the overall 
topological features remain unchanged. This reasoning 
explains that particle hole symmetry is relevant for the 
topology of the phase diagram only if the branching point 
is canceled by a zero of the on/off kinetics, i.e. precisely at 
K = K* . We conclude that the scenario for K = K* pre- 
sented in Fig. [23 generalizes the monomer TASEP cou- 
pled to Langmuir kinetics. Nevertheless there are many 
differences in the numerical values of several quantities, 
which are summarized in Table [I] 

An interesting difference between monomers and 
dimers resulting from the broken particle-hole symme- 
try can be read off from Fig. where the change of the 
phase diagram topology is shown as a function of ilr> 
and K . For monomers not only the critical value K* 
is shifted, but the line SId(K) is continuous. Both the 
lines £Id = and K — K* mark continuous transi- 
tions, but they concern two different events: while at 
K = K* the shape of the analytical solution changes, 
at SId = ftp a domain wall enters the system. In the 
monomers the critical detachment rates f2* and fll coin- 
cide and the topology of the phase diagram in Fig. [7\ is 
such that it is impossible to go from a configuration as in 
the bottom-left to the one in the top-right by changing 
K without crossing the transition line W D . The (contin- 
uous) transition occurring at firj = fl* D can be studied 
experimentally since the density profile are remarkably 
different. 

It is possible to generalize the continuous equation to 
the case of a system containing £-mers, using directly the 



continuous results presented in [l6L |2(| : 



p(l - tp) 



l-(£-l)p 



K 



(i-e P y 



(i-(i-i) p y- 



(36) 



The qualitative picture is contained already in the case 
of dimers; although the functional form of the solution 
depends on the size of the £-mers, the scenario of three 
branches, characterized by the maximum of the current- 
density relation p* and the horizontal asymptote pi, is 
robust. 

In conclusion we have studied a driven one dimensional 
lattice gas of dimers where the dynamics of the totally 
asymmetric exclusion process has been coupled to the 
on-off kinetics in the bulk. We used the ad-hoc "refined" 
mean-field for the TASEP part introduced in previous 
works 0, 0, 0| and proved that it is consistent with 
the assumptions made for the on-off part. The main ef- 
fect of extended nature of dimers on the phase behavior 
of the system is related to the breaking of a symmetry 
in the problem (particle- hole). This does have quanti- 
tative but not qualitative consequences on the density 
profile and on the phase diagram. The origin of the ro- 
bustness of the picture found for monomers can be traced 
back to the form of the stationary density profile which 
depends exclusively on the form of the current-density 
relation and of the isotherm of the on-off kinetics. This 
suggests that the TASEP dynamics washes out the inter- 
esting two-step relaxation dynamics that characterizes 
the on-off kinetics of dimers. The non-trivial outcome 
is that, in these systems, the diffusion (yet asymmetric) 
always dominates the large time-scale relaxation. Fur- 
ther studies on the dynamical correlation functions could 
point out more subtle differences between the dynamics 
of dimers and monomers. We conjecture that quite inter- 
esting phenomena would arise upon taking into account 
other interactions between the particles, which suggest to 
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Quantity 



Particle dens. 



Coverage dens. 



3 



Fa 



F d 



P(0) 



Jd 



PI 



K* 



Monomers 



O.Dp 



1-13 



a(l — a) 



0(1-/3) 



K 
1+K 



Dimers 



2p 



p(l-2p) 



QdP 



1-/3 
2 



(i+v^) 2 



a (1 — a) 
1 + a 



/3(l-/3) 
1+g 



1 



V1+4K 



TABLE I: Differences and analogies between the model with 
monomers and dimers: the current-density relation present 
in both cases a maximum at p* which appears at different 
values, but the particle-hole symmetry which is trivial in the 
monomers is actually effective in the dimers; the boundary 
conditions, as well as the current, need to be set differently 
in the dimers in order to keep track of their geometry; as a 
consequence of this the critical values of the boundary rates 
are different from the critical density in the case of dimers; 
the the isotherm is related to the binding constant in a non 
trivial way in the case of dimers; the particle-hole symmetry 
is preserved by the sources only in the monomers case. See 
the text for more details. 



can be expressed as 



p((*,0)|(*-l,2)V(t-l l 0)) = 

p((i-l,0)V(z-l,2)|(i,0))p(i,0) 
p((i-l,0)V(*-l,2)) 



(A2) 



The conditional probability in the numerator is unity, 
p((i — 1, 0) V (i — 1, 2)\(i, 0)) = 1, since a trail head cannot 
occupy site i — 1 if site i is empty. Moreover the events 
(i, 0) and (i, 2) in the denominator are mutually exclusive 
p((i-l,0)V(i-l,2)) = p(i-l,0)+p(i-l,2). Collecting 
results and using Eq. (|Alf> 



p(i,0|i- 1,0) = 



p(i - 1,0) +p(i- 1,2) 



(A3) 



The normalization condition compatible with the geom- 
etry of the dimers allows to rewrite the denominator, 
p(i- 1,0) +p(i- 1,2) = l-p(i-l,l) = l-p(i,2), and 
one obtains Eq. ©. 

A similar argument can be used to compute the joint 
probability p(i, 2; i+ 1, 0), using p(i + 1, 0|i, 2) under the 
same initial hypothesis 



p(* + l,0|*,2)=;p(t + l,0|(t,0)V(*,2)) 



(A4) 



Using this equation and the normalization condition p(i+ 
1, 0|t, 0)p(i, 0)+ P (i+l, 0\i, l)p(i, l)+p(i+l, 0\i, 2)p(i, 2) = 
p(i + 1, 0) (together with the fact that p(i + 1, 0\i, 1) = 
for obvious geometric reasons), one can show that p{i + 
l,0|i,2)~p(i+l,0|i,0). 



APPENDIX B: BOUNDARY CONDITIONS 



consider alternative current-density relations (with more 
maxima) or non-trivial on-off isotherm (cooperative on- 
off kinetics) . 
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APPENDIX A: THE CONDITIONAL 
PROBABILITY 

In this appendix we derive an approximate expression 
for the joint probability p(i, 0; i — 1, 0) and p(i, 2; i + 1, 0) 
starting from the mean-field assumption 

p(i, 0\i - 1, 0) ~ p(i, 0\(i - 1, 0) V (i - 1, 2)) . (Al) 

Since the event (i,0;i — 1,1) is excluded by the dimcr 
geometry, the mean-field approximation assumes the oc- 
cupancy of site i to be independent of the state of the pre- 
vious site. Using Bayes theorem, the r.h.s. of Eq. (|A1|) 



We recall the arguments given previously |3j, [ljj, [16| to 
justify the boundary conditions Eq. Ijl7|l . 

Left boundary: The left boundary conditions in the 
continuum limit can be derived from the rate equation. 
The incoming flux is given by the entrance rate multiplied 
by the probability of having the first two sites free, which 
is derived from the normalization condition on site i = 2 
(given that site i = 1 can never contain a trail head): 

p(l, 0; 2, 0) - 1 - p(2, 2) - p(2, 1) = 1 - p 2 - p 3 , (Bl) 

while the outgoing flux is determined by the usual rela- 
tion. This yields to: 



d T pi = a(l - p 2 - P3) 



P2(l -P3- Pa) 



1 



P3 



(B2) 



In the continuum limit (in first approximation) of the 
stationary state pi — pj, = p and therefore the equation 
can be solve for p in terms of a to get Eq. 117a|) 



P(0) 



1 + a 



(B3) 



This can be seen as if the left density was imposed by the 
reservoir density to which the system is coupled, which 
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is not p, but p/(l — p), therefore y- = j^— = a. Note 
that the current-density relation leads to the following 
boundary condition on the current: 



j(0)=3a 



a(l — a) 



1 



(B4) 



Right boundary: On the right boundary the last two 
sites are emptied at the same time and there is no ex- 
clusion on the last site i = N, therefore the incoming 
current on the second last site i = N — 2 is the usual one 
while the outgoing current is given by the occupation of 
the site i = N — 1 (since the last site is empty); consis- 
tently with the continuity property, the outgoing current 
at site i — N — 2 is also the gain term for the last site, 
while the loss term is given by fip^: 



drPN-l = 



Pn-2(1 - Pn-i - Pn) 



PN-l 



d T RN = Pn- 



Pn-i 



Pp 



N ■ 



(B5a) 



(B5b) 



It is obvious that in the continuum limit the density in 
the last two sites is different and only a coarse grained 
quantity like the average density (which is half the cov- 
erage density) makes sense: p = {pN-i + Pn)/2- One 



considers this last definition, the stationary state of 
Eq. HB5a[) with constant incoming current (that gives 
p(l — 2p)/(l — p) = pn-i) and the stationary state of 
Eq. IB 5 bp (that gives pn-i — Ppn), to obtain 



p(l - 2p) 2/3 



(B6) 



i-p r i-/r 

which finally gives the right boundary condition: 

p(l) = ^. (B7) 
The right boundary condition for the current reads: 



/3(l-/3) 
1 + /3 



(B8) 



which has naturally the same form of the left condition 
because of the effective particle hole symmetry. 

Note also that in the HD phase the particles wait a 
long time before detach giving rise to peculiar correla- 
tions: the (coverage) density profile exhibits a sawtooth 
profile superimposed to the analytical one. By using lin- 
ear analysis on the discrete map it can be shown that 
these correlations decay exponentially (but not as fast as 
the usual boundary layers) to the fix point. 
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